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Abstract. We study new FRW type cosmological models of modified gravity treated on the 
background of Palatini approach. These models are generalization of Einstein gravity by 
the presence of a scalar field non-minimally coupled to the curvature. The models employ 
Starobinsky's term in the Lagrangian and dust matter. Therefore, as a by-product, an ex- 
hausted cosmological analysis of general relativity amended by quadratic term is presented. 
We investigate dynamics of our models, confront them with the currently available astro- 
physical data as well as against ACDM model. We have used the dynamical system methods 
in order to investigate dynamics of the models. It reveals the presence of a final sudden 
singularity. Fitting free parameters we have demonstrated by statistical analysis that this 
class of models is in a very good agreement with the data (including CMB measurements) as 
well as with the standard ACDM model predictions. One has to use statefinder diagnostic 
in order to discriminate among them. Therefore Bayesian methods of model selection have 
been employed in order to indicate preferred model. Only in the light of CMB data the 
concordance model remains invincible. 
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1 Introduction 



As it is well-know a cosmological constant was the first and the simplest modification of 
general relativity performed by Einstein himself in order to stop cosmic expansion. Today 
the so-called Standard or Concordance (denoted also as ACDM or LCDM: Lambda Cold Dark 
Matter) Cosmological Model, which is based on this modification, turns out to be the best 
fitted model with respect to a huge amount of high precision currently available astrophysical 
data. In particular it properly describes a present day cosmic acceleration [1] by means of 
dark components: Dark Energy and Dark Matter. However ACDM model suffers for essential 
theoretical problems (e.g. well-known coincidence and fine tuning problems), specially related 
to a primordial stage of cosmic evolution [2]. Some of these problems can be cured by more 
sophisticated modifications. Among them the so-called /(i?)-gravity models constitute a huge 
family [3, 4] including also models based on Palatini formulation [5]. Quite recently there is 
a renewed interest in the Paltini modified gravity [6] which treats a metric and torsion-free 
connection as independent variables, see e.g. [5, 7, 8] for more details. Resulting equations 
of motion remain second order as in the Einstein gravity case. 

Modification of Einstein's General Relativity becomes a viable way to address the ac- 
celerated expansion of the Universe as well as dark matter and dark energy problems in 
modern cosmology (see e.g. [2] and references therein). This includes modified theories with 
a non-trivial gravitational coupling [9-13]. Viable non-minimal models unifying early-time 
inflation with late-time acceleration have been, particularly, discussed in [12]. Apart of cos- 
mological viability another justification for modified gravity should be taken seriously into 
account: curvature corrections appear naturally as a low energy limits of quantum gravity, 
quantization on curved background and/or effects from extra dimensional physics. 

Both purely metric as well as Palatini f{R)- gravity can be further extended by adding 
(scalar) field non-minimally coupled to the curvature [9-11]. In cosmological settings Palatini 
formalism gives rise, similarly to Einstein gravity, to the first order autonomous differential 
equation for a scale factor which can be recast into the form of Priedmann equation. This 
enables us to analyze the corresponding cosmological models as a one-dimensional particle- 
like Hamiltonian system with entire dynamics encoded in an explicitly obtained effective 
potential function. Such approach simplifies meaningly computer simulations and numerical 
analysis. Before doing these one has to constrain model parameters by using astrophysical 
data [14-17] and CosmoNest code [18]. Simple Palatini based cosmological models has been 
previously tested against various sets of cosmological data in [19]. Very recently it has 
been shown that gravitational redshift of galaxies in clusters is consistent both with general 
relativity and some models of /(i?)— modified gravity [20]. 

The paper is organized as follows. In Section 2 we recall some basic facts about Palatini 
gravity non-minimally coupled to dilaton-like field. Its cosmological application is considered 
in Section 3: particularly, general expression for generalized Friedmann equation is derived 
therein. Two classes of new cosmological models stemming from two different solutions of the 
so-called master equation are issued in Sections 4 and 5. In Section 6 we describe estimation 
of models parameters by astrophysical data followed by Bayesian method of model selec- 
tion provided in Section 7. Dynamical analysis of our models by means of one-dimensional 
particle-like Hamiltonian system with a Newtonian type effective potential function is sub- 
ject of our investigations in Section 8. In order to visualize the dynamics we watch plots of 
potential functions and draw two-dimensional phase space trajectories for the corresponding 
models. Section 9 treats about additional statefinder diagnostics. We end up with summary 
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of obtained results and general conclusions presented in Section 10. 

In this paper we took the first step in testing the kinematical sector of the cosmological 
models of modified gravity, i.e. we apply cosmographic analysis to the background dynamic. 
The next step will be the studying of matter perturbations in the class of models under 
consideration. In these future investigations the results of the present paper will be used as 
a starting point. This allows to exploit the advantages of Bayesian methods [21]. 

2 Preliminaries and notation 

The main object of our considerations in the present article is a cosmological application of 
some non-minimally coupled scalar-tensor Lagrangians of the type ^ 

L = ^{f{R) + F{R)Ld)+L^at (2.1) 

treated within the Palatini approach as in [9]. Hereafter we set = —i^g^'^ d^(j)du4> the 
Lagrangian for a free-scalar (massless) dilaton-like field (f) ^ and Lmat represents a matter 
Lagrangian. Because of Palatini formalism i? is a scalar R = R{g, T) = g^^'^ R^y{T) composed 
of the metric g and the Ricci tensor R^^iV) of the symmetric (= torsionless) connection 
r (for more details concerning the Palatini formalism see e.g. [7, 8]). Therefore (g', T) are 
dynamical variables. Particularly, the metric g is used for raising and lowering indices. 

We began with recalling some general formulae already developed in [9], both f{R) and 
F{R) are assumed to be analytical functions of R. Equations of motion for gravitational 
fields (r, (jr) can be recast [9] into the form of the generalized Einstein equations 

R^.u {h) = R^u (bg) = g^^P^ (2.2) 

(see also [8]), where R^y{h) is now the Ricci tensor of the new conformally related metric 
h = bg with the conformal factor b specified below. A (1, 1) tensor Pl^ is defined by 

and contains matter and dilaton stress-energy tensors: T,"?f* = ^^^siat ■ Tf., = Here one 

respectively has 

fc = (/ (i?) + F (R) La) = {L- L^at)/^ 
\b = f {R) + F' {R) Ld 

where a prime denotes the derivative with respect to R. 
Field equations for the scalar field <j) is 



(2.4) 



{^F{R)g^'''d^(t>) = (2.5) 

Dynamics of the system (2.1) is controlled by the so-called master equation 

bR = 2c- F{R)Ld + T (2.6) 

obtained by contraction of (2.2), where we set r = g^'^T^l^^ for a trace of the matter stress- 
energy tensor. In more explicit form it reads as 

2/ (R) -f'{R)R + T= {F' {R)R-F {R)) (2.7) 

These reproduce the same field equations as treated in [9] . 



^Throughout the paper we shall work with units c = SnG = 1. The metric signature is ( — h 
^One can easily add both mass and potential interaction for (f), cf. [9]. 
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3 Cosmology from the generalized Einstein equations 



Assuming the Cosmological Principle to hold we take the physical metric 5 to be a Friedmann- 
Robertson- Walker (FRW) metric g 

g = -dt^ + a^(t)\ dr^ + r^(de^ + sm^( e)d(p^]] (3.1) 

LI — Kr^ V /J 

where a = a{t) is a scale factor and k is the space curvature (k = 0, 1, —1). We also suppose 
that the matter content T^"* of the universe is described by a non-interacting mixture of 
perfect fluids. We denote by Wi the corresponding barotropic coefficients. Each species is 
represented by the stress-energy tensor T^*J = {pi +pi) u^Uy +Pig^v satisfying a metric (with 
the Christoffel connection of g) conservation equation V^^^^T^l] = (see [22]). This gives 
rise to the standard relations between pressure and density (equation of state) pi = WiPi and 

Let us recall that for the standard cosmological model based on the standard Einstein- 
Hilbert Lagrangian (considered both in the purely metric as well as Palatini formalisms) 

LeH = y/gR + Lmat (3.2) 

the corresponding Friedmann equation, ensuing from Einstein's field equations, takes the 
form 

+ ^.=\Y.^^-~'''^'"'' (3-3) 






when coupled to (non-interacting) multi-component non-interacting barotropic perfect fluids, 
where rjia~'^'^^~^^'^^ represents a perfect fluid component with an equation of state (EoS) pa- 
rameter Wi. Here H = ^ denotes a Hubble parameter. This is due to the fact that geometry 
contributes to the r.h.s. of the Friedmann equation through (cf. 2.7) 

R = -T = ^(1 - 3«;,)r?ia-3(i+-») = ^(1 - 3w^)pi (3.4) 

i i 

For example, the preferred ACDM model can be defined by three fluid components: the 
cosmological constant wa = —1, dust w^nst = and radiation u^rad = | assuming the 
spacial flatness condition k = 0. (As a matter of fact the spacial curvature term Ka~^ can 
be also mimicked by barotropic fluid Wcurv = ~3-) The radiation component which has 
no contribution to the trace r can be practically neglected due to extremely small value 
Clrad ~ 10~^. Alternatively, instead of introducing cosmological constant via perfect fluid 
Dark Energy component, one can, following Einstein, modify the Einstein-Hilbert Lagrangian 
(3.2): Leh Laeh = ^{R - 2A) + L^at- 

On the other hand we have that the field equation for the scalar field (p = cj){t) is 
■^{^F{R)4>) = 0, so that y/gF{R)4> = const and consequently gF{R)'^Ld = = const. 
This simply implies that 

F{R)^Ld = A^a~^ (3.5) 

with an arbitrary positive integration constant (see (2.5)). It means that this term behaves 
as a stiff matter component {w stiff = !)• 
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Assuming perfect fluid matter as a source, the generalized Einstein equations rewrites 
under the form 




^ F{R)Ld c {1 + 3wi)rji .(^^^^^ 

6b 6b ^ 6b ^ ' ' 

It becomes a generahzed Friedmann equation for the ordinary Hubble parameter H = ^ if 
we take into account that from (2.7) and (3.5) the scalars R and are implicit functions of 
the scale factor a. Further calculations give rise to the decomposition 

= K{a)G{a) (3.7) 

where 

G{a) = + J2 ii±|^a-^(i+-) - 2k(/' + F'L,)a- (3.8) 

and the function K{a) is defined by 



Kia) = 2{f + F'L, 



2/' - 4F'L, 



3[f + 2{F'R - F)Ld][f" + {F" - 2F~\F'Y)L, 



f'R. - /' + [F"R + 2F' - 2F-^{F'YR]Ld 

(3.9) 

where f = 'Ylii'^^i + 1)(1 ~ ■^Wi)rii a ^{i+'i'i) fQj. i^j^g multi-component fluid considered before: 
f = 3iff. As explained above the parameters, R, L^, and their functions, e.g. f{R), F'{R)L(i, 
etc. become through equations (2.7), (3.5) implicitly dependent of the scale factor a. In this 
way, the generalized Friedmann provides an autonomous system of first order ordinary differ- 
ential equation for the scale factor a (se the next Section for more details) . The decomposition 
(3.7) is furnished in such a way that for standard cosmology: b = l,c = R, = one has 
K{a) = 1 and one recovers the standard Friedmann equation (3.3), and particularly ACDM 
model as well. 

Our objective here is to investigate a possible cosmological applications and confront 
against astrophysical data of the following subclass of gravitational Lagrangians (2.1) 

L = ^(R + aR^ + PR^+^ + 7^'+"^d) + W (3.10) 



where a, /3, 7, 5, a are free parameters of the theory. It should to be observed that the grav- 
itational part f{R) contains the so-called Starobinsky quadratic term R? [23] with some 
R^^^ contribution. In the limit a , /? , 7 — )■ our Lagrangians reproduce General Relativ- 
ity. The constants a, (3, 7 are dimensionfull with the corresponding dimensions satisfying 
[R\ = [ai?^] = .. etc.. In fact, a numerical value of the constant 7 is unessential since it can 
be incorporated into the scalar field: — )■ |7|2i;^. Therefore, we further assume that it takes 
a discreet values 7 = 0, ±1. Such re-scaling would not be possible if one admits nontrivial 
self-interaction with nonvanishing potential energy f7((/>) for the scalar field. 

Following a common strategy particularly applicable within the Palatini formalism (see 
[7, 8, 10]) one firstly finds out an exact solution of the master equation (2.7). It allows to 
construct explicit cosmological model based on the generalized Friedmann equation. For this 
purpose and in order to reduce a number of independent parameters, we shall assume through 
the rest of this paper that universe is spatially fiat (k = 0) and filled exclusively with the 
most natural dust matter component only: f = p = — r = r]a~^ = 2/ — f'R + {F — F'R)Ld. 
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Modification comes from tlie geometric part of the theory. In this case formulae (3.7)-(3.9) 
simplify to more readable form 



2{f + F'Ld) [3/ - f'R + (3F - F'R)Ld] 



2f' 



^F'T. I ^\'if-f'R+{F'R-F)La][f"+(F"-2F~^{F')^)Li] 
/"R-/'+[F"_R+2F'-2F-i(F')^-R]id 



(3.11) 



This reconstructs the ACDM model under the choice f = R — 2A, F = 0, which is the limit 
a = 7 = 0, 5 = —1, /3 = — 2A. Einstein-DeSitter universe arising from the Einstein-Hilbert 
action (3.2) corresponds to a = /3 = 7 = 0. 



4 New cosmological models: solution I 

The first model we wish to be considered is based on the solution of (2.7), (3.5) resembling 
the Einstein-DeSitter universe (cf. 3.4) 

R = p = r]a~^ (4.1) 

provided that the integration constant A (see (3.5)) and parameter a take the values 



A' 



7 ? V , 



(4.2) 



where 6 7^ 0, 1. Positivity of can be ensured by an appropriate choice of 7 = ±1. 
Remember that for the ACDM case one has R = 4A + 770"^ instead (4.1). The conformal 
factor b reads now 



6=1 + 2ar]a~-^ + 



36-1 



-35 



(4.3) 



As a consequence, we have obtained the generalized Priedmann equations (3.7)-(3.9) 

= K{a)Gia) 



under the form 



and 



G(a) = -r?a-3 



3^ 



3<5 



,<5+l ^-3(5+1) 



(4.4) 



K{a) 



2 + Aaija 



,1-35 



-35 



2arja 



(l-35)(2-35) 



-35 



2 • 



1, 



5 

When a = and [3 = the Friedmann equation takes the simplest form H'^ = ^rya^", i.e. 
the same as in Einstein general relativity with a pure dust matter. As we already mentioned 
it reconstructs the (flat) Einstein-DeSitter model. 

Before proceeding further let us observe that the scaling properties of (3.6) are analogous 
to that in standard cosmology (3.3). It entitles us to introduce dimensionless cosmological 



(density like) parameters: the standard Oo,r 



as well the new one: 



0,/3 



^o,a = OiT] corresponding to new parameters a, (3 in the Lagrangian (3.10). With this set of 
variables the generalized Priedmann equation rewrites under the form 



\35 



2-2f]o,a(l + ^) = 



(l-3(5)(2-3(5) 
<5 



35 
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2 



(4.5) 



where redshift 1 + z = a ^. The number of free parameters (a, /3, S, r]) = (Oo,a> f^o,/3) ^, ^o,m) 
to be fitted by experimental data is 4. There are constrained by the following conditions: 
5/0,-1,1, f^o,me(o,i). 

Above equation can be rewritten in a more convenient for computer simulation form 



=no,mK{z)G{z) = no, 



2 + 40o,a(l + zf - 2i^l7o,/3(l + z) 



\3S 



2 - 217o,a(l + z) 
2{1 + zf + no,a(.l + zf 

Now normalization constraint can be simply set by 

Qo,mK{0)G{0) = 1. 



(l-3(5)(2~3(5) 



f]0,/3(l + 2) 



3<5 



3(5 



(4.6) 



(4.7) 



There is a number of interesting subcases, e.g, 5 = 4,1 which we shall not study here. How- 



ever the case of Einstein gravity supplemented by the quadratic Starobinsky term corresponds 
to /3 = 0: Les = \/g{R + aR'^) + imat- In this simplest case one gets 



\^oJ [2 - 2J7o,a(l + z)^] 



(4.8) 



Cosmology for (4.8) has been also recently investigated in [24] (see also [25]). We will see by 
comparison against experimental data that this model does not differ qualitatively from the 
ACDM, i.e. as expected it contains dynamical cosmological constant which becomes active 
in recent times. Adding dilaton field (j) makes it even more quantitatively similar provided 
we shall use the same standard solution (4.1). 



5 New cosmological models: solution II 

Yet another models can be determined by the solution 



1] 



(1-6) (3 



a = 25 



(5.1) 



where ^ = ^^^^ 
The conformal factor b reads 



1 

i+s 



One needs /3(1 — S), ^5 > 0. Now the parameter /3 cannot vanish. 
1 + 4(5 



26 



+ 2a^a~~ + /3(1 + 6)Ca~~ 



(5.2) 



As a consequence, we obtain components constituting the generalized Friedmann equations 
(3.7)-(3.9) (cf. (3.11)) under the form 



G(a) 



1 + 6 



1 .2 



2-6 



36 ^ 3 ^ 3(1 - (5) 



77a 
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1+4(5 



Kia) 



+ Aa^a'T+s + 2/3(1 + 6)^^a~~^ 



1+4(5 



_36_ 
+S 



5 ^ ^ 1+5 

Introducing, as before, dimensionless density parameters: i^o,m 



(5.3) 



3m 



3m 



(l-(5)/3 



1 



0,Qf 



aH^Qo^ls one gets in terms of the redshift 1 + z = a ^ 



i7 



K(z)G(z) = 
'1 + 5 



1+4(5 



O0,/3(l + Z) 1+^ + 30o,a!^0,/3(l + Z) i+« + ^— |Oo,m(l + zf 



(5.4) 



The number of free parameters (a, /?, 6, rj) = (fio,a) ^o,/3) 5, ^o,m) to be fitted by experimental 
data is 4. There are constrained by the following conditions: (1 — (5)/3, r/, r^o,/? > 0, ilo,m £ 
(0,1), 5/0,-1,1. 

Introducing f]o,c = ^^o.mf^o ^ allows us to rewrite the last equation in a more convenient 

form 



(- 



1+4(5 
(5 



+ 12J]o,a(l + z)^ + 2i±f Oo,c(l + z) 



1+S 



1 



1+4(5 



3i 



+ 6f^0o,„(l + z)^ + fE|^^o,c(l + z)-^ 

2-5 



-(1 + z)— + 30o,a(l + 2)— + Y37^o,c(l + ^)' ^5.5) 



Then normalization constraint can be simply set by 

J7o,/3i^(0)G(0) = 1. 



(5.6) 



Again in order to better control the role of Starobinsky term one can switch it off by setting 
q; = 0. As it has been mentioned before the parameter (3 cannot vanish now. 

One can summarize this part by concluding that we have obtained for numerical analysis 
four new cosmological models which will be further on denoted correspondingly as /, Ia=o 
and II,IIa=Q. Models / and // have three free parameters to be fitted by experimental 
data. Models Ia=o IIa=o have two such parameters. These models are formulated not 
only in terms of Lagrangian functions but also by giving explicit form for the corresponding 
Friedmann equations. Finally, the simplest case //3=o describes cosmological model based on 
Palatini modified gravity without dilaton field which has similarly to ACDM only one free 
parameter. Such quadratic gravity models have been extensively studded in the literature 
(4.8) by using different methods than the one we have employed here. The reason is that the 
explicit form of the corresponding Friedmann equation for this case (see (4.8)) has not been 
used except [24]. 



6 Constraining model parameters by astrophysical data 

In the Bayesian approach to estimation model parameters (i.e. best fit values and credible 
intervals) one uses a posterior probability density function (pdf), which is defined in the 
following way: 
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Q is the vector of model M parameters and P{D\Q, M) = L is the hkehhood function for 
model M. The so called prior pdf for model parameters P{Q\M) should be estimated before 
the data D comes into analysis and should involve information which we have gathered in 
earlier studies, e.g. with different data sets or on theoretical grounds. Finally P(D\M) = E 
is the so-called evidence, which could be ignored in model constraint analysis. The best fit 
values can be estimated using the maximum of the joined posterior pdf (6.1) (i.e. its mode). 
One can also consider marginalized posterior pdf: 



p{e,\D,M) = / p{e\D,M)d4>, 



(6.2) 



where 6> = {^j, c^}, and use its mode or mean and credible interval (usually defined as interval 
which involves 68% or 95% of the probability) . 

In order to estimate the parameters of our models we use supernovae (SNIa) data [14], 
the observational H[z) data [15], the measurement of the baryon acoustic oscillations (BAO) 
from the SDSS luminous red galaxies [16] and information coming from CMB [17]. 

We use a sample of = 557 data [14], which consist of data from Union [26], Supernova 
Cosmology Project [14], SDSS SN Survey [27] and CfA3 [28]. The hkelihood function is 
defined in the following way: 



LsN oc exp 



theor i,obs"i 



2af 



(6.3) 



where: fjj is the total measurement error, ^f^^ = rrii — M is the measured value (mj-apparent 
magnitude, M-absolute mag nitude of SNIa), = 5log-^o Dli + M = blog^Q dLi + 25, 

Ai = — 51ogioffo + 25 and Du = Hodu, where du is the luminosity distance given by 
dii = (1 + Zi)c /q^' h(z') (witli the assumption k = 0). In this paper the likelihood as a 
function independent of Hq has been used (which is obtained after analytical marginalization 
of formula (6.3) over -f^o)- 

We use constraints coming from 13 measurements of Hubble function at different red- 
shifts z. Those data points are obtained by two observational methods. The first method is 
based on the measurements of spectroscopic ages of red galaxies [29] , while the second one 
is based on the measurements of BAO scale in radial direction [30]. For the H{z) data the 
likelihood function is given by: 



Lh^ oc exp 



{Hizi 



2aj 



where H{zi) is the Hubble function. Hi denotes observational data. 

We also use information coming from the so called BAO A parameter, which is related 
to the Baryon Acoustic Oscillations scale measured in the redshift space power spectrum of 
luminous red galaxies (LRG) from the Sloan Digital Sky Survey (SDSS) [16]. For BAO A 
parameter data the likelihood function is characterized by: 



Lb AO oc exp 



^^thcor y^obs^2' 



(6.4) 



where = 



m,0 



H(za) 
Ho 



ZA JO 



2A Hp 



dz 



and A"^"" = 0.469 ± 0.017 for za = 0.35. 



-9- 







S 


^0,/3 


no.™ 


Xtot/2 




= 0, n 


0,c 


e< - 


1, 10 >, S £ (0, 1 > 


2.572 


0.997 


0.254 


0.652 


294.233 




= -100, 


no 


c 6< 


-1, 10 >, i5 6 (0, 1 > 


1.748 


0.553 


0.012 


0.045 


278.974 




= -150, 


no 


c e< 


-1, 10 >, <5 e (0, 1 > 


2.226 


0.538 


0.005 


0.024 


278.853 




= -300, 


no 


c e< 


-1, 10 >, <5 e (0, 1 > 


3.469 


0.524 


0.009 


0.013 


278.737 


Ho. a 


= -1000, 


no 


,c e< 


-1, 10 >, (5 s (0, 1 > 


4.693 


0.508 


0.020 


0.019 


278.655 



Table 1. Comparison of estimated parameters r2o,c, ^^o,/3 and f2o,m for model II with fixed value 
of rio,Q- It shows that essential parameters as S or rio.m behaves stable under a wide range of ilo.a 
provided IIo.q 7^ 0. 



Finally, we use constraints coming from CMB temperature power spectrum, ie. CMB 
R shift parameter [31], which is related to the angular diameter distance {Da{z^,)) to the last 
scattering surface: 

j^^Vn^ (6.5) 

c 

The likelihood function has the following form: 

I {R — RohsY 



LcMB oc exp 



2 a\ 



(6.6) 



where Robs = 1-725 and aj^^ = 6825.27 for = 1091.3 [17]. 
The entire Ltot is characterized by: 

Ltot = LsnLh,LbaoLcmb- (6.7) 

We assume flat prior probabilities for model parameters (second column of Table 3) . The 
prior probabilities because of constraints coming from previous estimations of the parameters 
were calculated using N = 192 SNIa data sample [32]. Additionally we have assumed that 
Hq = 74.2 [kms~^Mpc~^] [33] and 0^ G [0, 1]. 

The knowledge of prior pdf for model parameters and the formula for likelihood function 
enable us to calculate joined and marginalized posterior pdf. The mode of joined posterior 
pdf as well as mean (together with 68% credible interval) of marginalized posterior pdf were 
calculated using publicly available CosmoNest package [18, 34-36], which was modified for 
our purpose. 

In table 3 we display the best fitted parameters for all our models Iq,=o, -f/3=0) IIa=o, II 
as well as for ACDM estimated by CosmoNest package. Parameters r^c^S) ^o,m for models I 
and parameter Oo,™ for models II are calculated using equations (4.7), (5.6) and estimated 
values of the remaining parameters. We consider two cases: estimations with data sets com- 
ing from late universe (i.e. SNIa, H(z) and BAO) and estimations including also information 
from early universe (i.e. SNIa, H(z), BAO and CMB). Top part of the table relates to the 
first case, while the bottom relates to the second one. The best fit values correspond to the 
mean of the marginalized posterior pdf (68% credible intervals are also shown, estimation of 
sample variance is shown in cases where the credible interval has misleading values). Finally 
mode of joined posterior pdf is given (values in brackets) and the corresponding value of 
(table 2) (see next Section for explanations). 

In order to get better insight into the role of quadratic term we have made several 
estimations for fixed values of Oq. The results are gathered in table 1. They differ significantly 
when a = and a 7^ 0. In the last case the concrete values of Vta have secondary meaning. 
Moreover, the switching off R? term leads to the cosmology which is not longer similar to 
ACDM (cf. Fig 11, 9). 
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Figure 1. Constraints of tiie parameters of model / - estimations without CMB data. In 2D plots 
solid lines are the 68% and 95% confidence intervals from the marginalized probabilities. The colors 
reflect the mean likelihood of the sample. In ID plots solid lines show marginalized probabilities of 
the sample, dotted lines are mean likelihood. For numerical results sec Tabic 3, No 3. 



At the fig. 1, 3, 2, 4 we show constraints of the parameters of our models. In the 2D plots 
solid lines are the 68% and 95% confidence intervals from the marginalized probabilities. The 
colors reflect the mean likelihood of the sample. In the ID plots solid lines show marginalized 
probabilities of the sample, dotted lines are mean likelihood. It is seen that it's not possible 
to constrain properly parameter r2o,a using available data (cf. table 1). 

It is interesting to compare (Fig 5) a posterior probability density functions of Oo,m 
for new models with the one for ACDM model (see table 3 for the corresponding best fit 
values and credible intervals). As one can conclude only for model / (all three parameters 
fitted with and without CMB data) it remains at the same level when compared with ACDM 
model. It should be remarked that in some cases, best fit value for Oo,m is small and close 
to the well-known amount of barionic matter 0^ ~ 0,05 (cf. table 3 row no 1, 5, 10, 11). 
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Figure 2. Constraints of the parameters of model / - estimations including CMB data. In 2D plots 
solid lines are the 68% and 95% confidence intervals from the marginalized probabilities. The colors 
reflect the mean likelihood of the sample. In ID plots solid lines show marginalized probabilities of 
the sample, dotted lines are mean likelihood. For numerical results see Table 3, No 9. 



Finally, at the fig 6 Bubble's diagrams for all our models are shown and compared 
with ACDM. Most of them are practically indistinguishable from the concordance model on 
the level of Hubble diagrams. Moreover, the x^~test gives mostly comparable results and 
in few cases priority belongs to new models provided that CMB data are neglected. It is 
clearly seen that only for the special case IIa=o (fitting made with and without CMB data), 
the corresponding plots considerable differ from all remaining ones. These motivates us for 
further analysis and employing Bayesian model selection methods. 

7 Models comparison and selection 

Bayesian theory enables to compare investigated models, i.e. enables to show which one is the 
best (most probable) in the light of analysed data. Let us consider the posterior probability 
for model indexed by i (Mi): 

P(M.|«).i»M. (7.1, 
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Figure 3. Constraints of the parameters of model // - estimations without CMB data. In 2D plots 
solid lines are the 68% and 95% confidence intervals from the marginalized probabilities. The colors 
reflect the mean likelihood of the sample. In ID plots solid lines show marginalized probabilities of 
the sample, dotted lines are mean likelihood. For numerical results see Table 3, No 5. 



P{Mi) is the prior probability for the model under investigation, P{D) is normalization 
constant, P{D\Mi) is the marginalized likelihood (also called evidence) and is given by: 

P{D\Mi) = j P{D\Q,Mi)P{Q\Mi)dQ. (7.2) 

It is convenient to choose the base model and compare all models from investigated 
set with respect to this base model by considering the ratio of posterior probabilities (posterior 
odds): 

P{Mb\D) ^ P{M,) P{D\Ah) ^ PjM,) 
P{Mi\D) P{Mi)P{D\Mi) P{Mi) 

When all considered models are equally prior probable the posterior odds is reduced to 
the evidence ratio, so called Bayes Factor {Bbi). Its value can be simply interpreted as 
the strength of evidence in favour of base model : < InS^j < 1 as inconclusive, 1 < 
lni?5j < 2.5 as weak, 2.5 < InSfoj < 5 as moderate and InSbj > 5 as strong evidence. Since 
In i?5j = — In Bi^i the negative values of In Bi^i should be interpreted in favour of model under 
investigation. 
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Figure 4. Constraints of the parameters of model // - estimations including CMB data. In 2D plots 
solid lines are the 68% and 95% confidence intervals from the marginalized probabilities. The colors 
reflect the mean likelihood of the sample. In ID plots solid lines show marginalized probabilities of 
the sample, dotted lines are mean likelihood. For numerical results see Table 3, No 11. 



In Table 2 one can find values of logarithm of Bayesian Factor for new models, which 
were calculated with respect to the base ACDM model using CosmoNest code. The values 
were averaged from five runs. 

As one can generally conclude models based on solution II are not supported by data 
used in analysis: we found weak and strong evidence in favor of ACDM model (for both with 
and without CMB data cases). When considering Model / one can find out that there is 
weak evidence in favor of it with respect to ACDM model in the light of observations from 
late Universe. When information from early Universe is included, i.e. CMB data point is 
taken into account, the conclusion change: moderate evidence is found in favor of ACDM 
model. Amazingly, the special case of quadratic gravity (/^=o) is not supported by the data, 
in both cases the evidence in favor of the base model is strong. The result for model Ia=o is 
inconclusive when the data from late Universe are considered. Admitting information from 
CMB gives us strong evidence in favor of ACDM model. 
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Figure 5. Posterior probability density functions of ftm,o parameter for all cases (red lines). First 
row correspond to models I. second to models II. First column and forth: P ~ (i.e. quadratic 
gravity); second and fifth: a — . In third and sixth column 3 parameters were fitted. Black curves 
correspond to ACDM model. Left panel: parameters were fitted using Snia, and BAO data. Right 
panel includes CMB data. In some cases ri„i,o is the same order as Qbarionic ^ 0.04. For numerical 
values see Table 3. 




Figure 6. Comparison of Bubble's diagrams for model I (green) and II (magenta). Grey line 
denotes special case of quadratic gravity Ip=o. Blue (/q,=o) a-nd light blue (//q=o) lines denote most 
divergent with respect to ACDM (black) models without Starobinsky's term. Left panel corresponds 
to estimation without CMB data, right panel relates to estimation including CMB data. 

8 Phase space descriptions of model dynamics 

Having fixed free parameters of our models on can think about their dynamical properties 
encoded in a time evolution of the scale factor. In fact, the dynamics of our models is deter- 
mined by very complicated Friedmann type equations which constitute first order ordinary 
(non-linear) autonomous differential equations on the scale factor. Moreover r.h.s. of the 
generalized Friedmann equation (4.5,5.4) is rather rational than polynomial function (cf. 
(3.3)) of a. While obtaining exact solutions of these equations is very difficult (if possible 
at all) it is enough to apply qualitative and numerical methods of analysis of differential 
equations. The main goal of this method is instead of studying individual trajectories of the 
system under consideration to describe a geometrical structure of a phase space. The main 
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Model 


Estimation without CMB data 


Estimation with CMB data 




In B^c DM, Model 




In B\c DM, Model 




la=0 


-0.8 ±0.4 


272.287 


453.3 ±0.6 


706.235 


1/3=0 


1665.2 ±0.3 


274.142 


35.5 ±0.2 


310.949 


I 


-2.2 ±0.3 


271.400 


4.2 ±0.3 


276.668 


IIa=0 


21.3 ±0.2 


294.233 


169.1 ±0.3 


438.839 


II 


6.5 ±0.3 


279.237 


206.6 ±0.3 


475.881 


ACDM 





276.583 





276.726 



Table 2. Values of the logarithm of Baycsian Factor together with the corresponding for models 
based on solutions I and II, with respect to ACDM model. 



advantage is the possibility of investigating entire evolution, represented by trajectories in 
the phase space, for all admissible initial conditions. The phase space (the phase plane (a, d) 
in our case) is organized through critical points and trajectories. As a result we obtain a 
global phase portrait of the system illustrating the stability of special solutions as well as its 
generic properties among all evolutional paths in the phase space. 

It would be useful to represent the dynamics of the cosmological models in terms of 
dynamical system theory. For this aim let us re-parametrize the original time variable in 
(4.5,5.4) to the new re-scaled variable, say r such that dr = li^olc^^ and define an effective 
potential function V{a) = —^a^K{a)G{a). Then we obtain that dynamics is reduced to the 
dynamics of a fictitious particle of unit mass moving on a half line in the effective potential 
V{a) with energy level = —k/2, where k = 0, ±1 is the spatial curvature index (for 
generality one considers models which are not flat): 



1 f da^^ 



2\dT 



+ V{a) = (8.1) 



The system is defined in the part of the configuration space in which E^^ — Vis non-negative. 
Therefore the dynamics of cosmological model is governed by dynamical system 

Tr-^ (8.) 
t-fa 

of Newtonian type. Critical points of this system correspond to extremes on a diagram of the 
effective potential. Because the stability of the critical point depends on the eigenvalues of 
a linearization matrix (solutions of the characteristic equation) calculated at this point one 
obtains two possibilities (a trace of the linearization matrix is always vanishing). If one has 
a maximum then the critical point is a saddle. In contrast a minimum on the diagram of the 
potential function means that we have a center type of critical point. Critical points provide 
stationary solutions: stable for minimum and unstable otherwise. 

In fact dynamics of the system can be read off from the graph of potential function itself 
(in cosmological setting the effective potential is non-positive: 1^ < 0). Since conservation 
of (total) energy holds, each trajectory is labeled by some energy level E = constant (might 
be negative). The difference E — V '>Q \s equal to kinetic energy and provides information 
about the velocity of cosmic expansion: increasing values of the potential are slowing down 
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the evolution. Because of this the region below potential plot is forbidden for the motion. 
Particularly, points on the graph are turning points. For example, models with a classical 
bounce contract to some minimal value of the scale factor and then after approaching this 
value expand (to infinity). Oscillating models are, therefore, result of a potential cavity. As 
a result one can simply obtain different types of evolutions admissible by different values of 
energy levels. For example, as it was mentioned above, different choice of spatial curvature 
index k will determine different cosmic evolution. Such visualization is complementary to 
the one which offers phase portrait. 
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Figure 7. The diagram of the effective potential in particle-like representation of cosmic dynamics 
for model I versus ACDM model (left picture relates to estimation without CMB data, the right 
relates to estimation employing CMB data; table 3, No 3, 6, 9, 12). Note that till the present epoch 
two potential plots almost coincide. Particulary, one can observe decelerating BE era. Maximum of 
the potential function corresponds to Einstein's unstable static solution (saddle point). Discrepancies 
become important in the future time: e.g. discontinuities of the potential functions (vertical, red 
lines) denote that V — > — oo, i.e. a — >■ oo for a — >■ a-''™"'. It turns out to be finite-time (sudden) 
singularity, hi any case the shadowed region below the graph is forbidden for the motion. 
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Figure 8. The diagram of the effective potential in particle like representation of cosmic dynamic for 
the model of quadratic gravity /^=o versus ACDM model (left picture relates to estimations without 
CMB data, the right relates to estimations employing CMB data; Table 3, no 2, 6, 8, 12). Maximum of 
the potential function corresponds to unstable static solution (saddle point). Again, until the present 
epoch there is no striking differences between plots. One can observe finite-size sudden singularity in 
the near future (vertical, red lines). In any case the shadowed region below the potential is forbidden 
for the motion. 
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Figure 9. The diagram of the effective potential in particle like representation of cosmic dynamic 
for the model Ia=o versus ACDM model (left picture relates to estimations without CMB data, the 
right relates to estimations employing CMB data; Table 3, no 1, 6, 7, 12). Both pictures represent 
decelerating, till the present epoch, universe with BB initial singularity. Now discrepancies with 
ACDM model are more evident. In any case the shadowed region below the potential is forbidden for 
the motion. 
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Figure 10. The diagram of the effective potential in particle-like representation of cosmic dynamics 
for model II versus ACDM model (left picture relates to estimation without CMB data, the right 
relates to estimation with CMB data; Table 3, no 5, 6, 11, 12). The evolution of the model is 
represented through the energy level. Therefore both models arc bouncing type. The universe is 
contracting, reaches the minimal size and then is expanding with acceleration. Notice that some part 
of the potential plot coincide with ACDM model. 




Figure 11. The diagram of the effective potential in particle like representation of cosmic dynamics 
for the model IIa=o versus ACDM model (left picture relates to estimation without CMB data, the 
right relates to estimation with CMB data; Table 3, no 4, 6, 10, 12). One can observe big discrepancies 
with respect to the concordance model. Moreover, for the first time, adding CMB data changes the 
type of expansion from accelerating to not accelerating. 
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On the diagrams we have illustrated the effective potential functions (Fig. 7, 9, 8, 10, 
11) as well as generic samples of phase portraits for our models (Fig. 12, 13). We found 
that although analytical formulae for effective potential functions for various model are very 
different their numerical values and therefore a shape of the corresponding graph might be 
almost the same. As it was explained above similar shapes give rise to similar dynamical 
behavior. This is what happens in our case. 

V[a] 




Figure 12. The diagram of the effective potential and the corresponding phase portrait of the model 
/ for estimations using CMB data. However for estimations without CMB data the phase portrait 
looks similar. We marked as a bolded trajectory of the flat model determined by the energy constraint 
E = 0. The phase space is divided by this trajectory on two domains at which lie closed and open 
models. We have situated the only possible critical point, the saddle (astatic, 0) which represents the 
Einstein static universe. The vertical red line a = astatic passing through the saddle critical point 
divides each trajectory into two parts: decelerating (V{a) is a growing function of its argument) and 
accelerating {V{a) is a decreasing function of the scale factor) eras. This part of the phase portrait 
is topologically equivalent to phase portrait of ACDM model. Particulary, the Bing Bang era is 
decelerating. 
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Figure 13. The diagram of effective potential and the corresponding phase portrait on the plain (a, a) 
for model // - estimations including CMB data. (In the case without CMB data the phase portrait 
looks similar.) The phase portrait shows all evolutional paths for all admissible initial conditions. We 
mark as bolded the trajectory of the flat model determined by the energy constraint E = 0. The 
phase space is divided by this trajectory on two domains at which lie closed and open models. Note 
that all solutions are bouncing type. In contrast to ACDM model on can notice accelerating scenario 
from the very beginning. 



Let us discuss the plots of the effective potentials for new models in more details. The 
model / has reached very good numerical agreement with ACDM for a wide range of the scale 
factor a: from the very beginning till the present time (fig. 7). After Big Bang the universe 
expansion slows down (deceleration epoch) till some critical value (equal to static solution) 
and since then speeds up. This cosmic acceleration lasts forever for the case of ACDM. 
In contrast, the model / predicts finite size future singularity since the effective potential 
function (as well as velocity plot) has a pole at a^j"'°'^ = 2.359. It can be shown that this pole 
has to be reached at finite time. Therefore it is known as a sudden future singularity [38]. 



- 20 - 



In the case of quadratic gravity (/^=o)) the shape of the potential function is the same as for 
model / with poorer numerical coincidence with ACDM (fig. 8): of*""' = 1-639. The initial 
deceleration era characteristic for ACDM model is absent only for the case of model // which 
offers a model of permanently accelerating universe of bouncing type. After a short period of 
strong acceleration (possible inflation) there is an era of low acceleration followed by ACDM 
epoch. A characteristic for models / and ACDM saddle point is not present for the case of 
model //. Finite size future singularity of the potential function appears at Oj™"' = 2.343. 
However on the Hubble diagrams differences between ACDM and models /, //3=o and II, are 
negligible. Similarity to ACDM is lost for both models with a = (see fig. 9, 11) which 
makes this case difficult to accept. Models of modified gravity enlarge zoo of possible cosmic 
evolutions because of the presence of finite time final size singularity of the potential function 
at some finite value of the scale factor. 



9 Other analysis and diagnostics 

In this paper we have studied and confronted against astrophysical data as well as against the 
standard ACDM model five cosmological models obtained from two solutions presented in 
Sections 4, 5. For a deeper examination and comparison of properties of our models we have 
used various cosmological parameters as: deceleration parameter q{a), effective equation of 
state Weffia), JERK j{a), SNAP s(a) (see e.g. [39, 40]). For a later convenience we shall 
express them in terms of the effective potential V{a) 
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where V'{a) = Therefore, explicit analytic form of the potential function is also helpful 
for determination of these diagnostics. Below we shall summarize the results of our analysis, 
for details see fig. 14, 15, 16. 

Deceleration parameter q{a) for ACDM, in contrast to other models, goes asymptotically 
to —1. However, the present day values of q{a) for new models are not much different, 
particulary qi^^^^^a = 1) < qi{a = 1) < qACDM{cL = 1) < qn{a = 1). Only for two models 
with a = the present day values o q{a = 1) are visibly greater (see fig. 14 and table 3). 
One notices that the values of q{a) are almost the same at the point a ~ 0.75 for all of the 
models except the ones with a = 0. At the beginning q{a = 0) ~ 0.5 for all models except 
// for which one has qn{0) = —0.02 instead. It means a matter dominated epoch at the 
beginning of cosmic evolution (cf. below for plots of Weff{a)). 

Very similar situation is for effective equation of state parameter since it is algebraically 
related to q{a) (cf. (9.2)). For ACDM model lima^-oo it'e//(«) = —1 and Weffj^^gia = 1) < 

Weff,l{a = 1) < Weff,ACDM{a = 1) < W^ffjlia = 1) < U)e//,7„=o(« = 1) < 'f«e//,//a=o (" = 1) 
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Figure 14. Plots of the deceleration parameter q{a) (left panel) and the effective equation of state 
Wef f (a) (right panel) for all models under investigation. Only model // (magenta) provides permanent 
acceleration. Models Ia=o (blue) and IIa=o (hght blue) have no acceleration epoch at all. There is 
intriguing intersection near a = 0.75 for plots representing four models: /, //3=o, II and ACDM. Thin, 
vertical line denotes present time. The model parameters were fitted using SN, Hz, BAO and CMB 
data (Table 3, No 7-12). 




Figure 15. Plots of parameters JERK j{a) (left panel) and SNAP s{a) (right panel) for investigated 
models. From the figure one can see that different models predict different prcscnt-day values of 
JERK and SNAP. Unfortunately estimations of these parameters arc beyond our present observational 
possibilities. However some recent analysis support the current values of jerk bigger than 2 [40]. Thin, 
vertical line denotes present epoch. The model parameters were fitted using SNIa, , BAO and CMB 
data (Table 3, No 7-12). 

(see bottom fig. 14). At a = there is Weff,ip^Q ~ ^e//,7 ~ Weff,ic,=o ~ Weff, 11^,^0 ~ 
Wef f,\CDM = 0. This explicitly indicates early matter dominated epoch noticed already 
before for these models. 

As it can be seen on the fig. 15, that JERK parameter j(a) for ACDM model is 
constant Jacdm = 1 and for a < 0.5 models /q=o, I and I 0=0 approximate this value. 
Similarly, jua^o ~ 0-2 is also constant. From the figure one can see that different models 
predict different present-day values of JERK and SNAP. Unfortunately estimations of these 
parameters are beyond our present observational possibilities. However some recent analysis 
support the current values of jerk bigger than 2 [40]. 
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Figure 16. Plots of parameter Om{z) allowing to answer the question on nature of dark energy 
[41]. For ACDM (black line) we have Om{z) = constant which means that cosmological constant 
is the best description of dark energy. For models /a=o (gray) and /q=o (blue) the function Om is 
increasing what means that the best model for dark energy is provided by a phantom. Remaining 
models have decreasing Om functions and therefore they are of quintessence type. However, one can 
notice that both / and /^=o models have a long period of being of cosmological constant type. The 
model parameters were fitted using SN, H^, BAO and CMB data (Table 3, No 7-12). 

In order to better understand properties of our models we shall use Om diagnostic from 



which is designed to answer the question of nature of dark energy. For ACDM model we have 
Om{z) = constant (see black line on fig. 16) which means that cosmological constant A is 
the best description of dark energy. For the model 1/3=0 the function Om is increasing what 
means that the best model for dark energy is provided by a phantom. Remaining models 
have decreasing Om function and therefore they are of quintessence type. However, one can 
notice that two models / and 1/3=0 have a long period of being of cosmological constant type. 

10 Summary and conclusions 

In this paper we have found, and then analyzed by fitting to experimental data, two families 
of cosmological models based on two different solutions of Palatini modified gravity equipped 
with non- minimal curvature coupling to (free) scalar dilaton-like field. We assume Cosmolog- 
ical Principle to hold and standard spatially flat FRW metric with dust matter as a source. 
Our analysis reveals how rich is class of cosmological models offered by modified gravity. Be- 
side of two main models denoted respectively as / and // we have investigated three special 
cases 1/3=0, Ia=o and IIa=o corresponding to reduced lagrangian functions. Particularly, as 
a by-product, we have employed quadratic a la' Starobinsky gravity, I/3=o, for description 
of the cosmic acceleration. All models labeled by / correspond to the solution (4.1) which 
is the same solution as in Einstein gravity provided with FRW metric and dust matter. In 
the case of standard gravity this solution leads, after solving Friedmann equation, to decel- 
erating expansion. Modification of gravitational Lagrangian (e.g. by adding cosmological 
constant) modifies Friedmann equation and thence give rise to accelerating expansion. In 
contrast, models labeled by // are based on the solution (5.1) which has no correspondence 
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z>0 . 



(9.5) 
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to standard gravity. Therefore the hmit /3 i— )■ cannot be performed. Having estimated 
models parameters by using recent astrophysical data one has demonstrated that different 
models become in a good agreement with the data. We have also performed various com- 
parative analysis against ACDM concordance model which is commonly accepted as a best 
fitted model. 

The first model /, in spite of much more complicated Friedmann equation, qualitatively 
and quantitatively mimics the ACDM one (f^oA = 0.77) from the very beginning of cosmic 
evolution (Big Bang) until the recent time. As it can be seen, on various plots, effective 
potential, deceleration and effective equation of state parameters (fig. 7, 14) are almost the 
same until the present time (a = 1). Differences appear for diagnostics employing third and 
forth order derivatives of the scale factor (fig 15). These properties would be helpful for future 
discrimination between models. Both models exhibit existence of the initial singularity, but in 
contrast to ACDM, our model predicts the (final) finite-size sudden (finite-time) singularity 
at the point a = 2.357. The effective potential for our model has one maximum at the point 
a = 0.552, {z = 0.811) which corresponds to saddle critical point at the phase portrait (see 
fig 12) and representing Einstein's static solution in full analogy with ACDM. For the present 
time the deceleration parameter q = —0.814 and WQ^eff = —0.876 are within the expected 
estimations. The Hubble diagram for this case is practically indistinguishable from ACDM 
diagram, cf. Fig. 6. Moreover, the estimation with CMB data added does not provide any 
essential (qualitative nor quantitative) changes to this model. 

At the fig. 5 are shown posterior probability density functions for f^o.m parameter for 
all of our models. It's clearly seen that density function for model I almost entirely overlaps 
such function for ACDM model (see also table 3 for the corresponding best fit values of 
^o,m)- One can also see that inclusion to our estimations CMB data makes possible proper 
estimation of Oo^m for the following cases: I/3=o, IIa=o and II. In particular, for models II 
amount of Oo.m is one order less than for ACDM and it has the order of barionic matter. 

The second our model gives rise to permanently accelerating Universe without initial 
singularities: Big Bang scenario is replaced by Big Bounce. The potential plot has no max- 
imum and there is no saddle point at the phase portrait (fig 12). The lack of the initial 
singularity is the main advantage of this model in comparison with ACDM model. The coin- 
cidence with ACDM effective potential holds for a € (0.6, 1.2). The acceleration parameter 
is almost linear for the scale factor near zero. The potential function indicates big slope for 
the scale factor near zero what can be interpreted as inflationary era. Effective equation of 
state at the beginning is close to — | (fig 14), which corresponds to domination of the spatial 
curvature at the initial phase of Universe's evolution. For the present time the deceleration 
parameter qq = —0.529 and effective equation of state wo^eff = —0.686. Generic finite time 
finite size singularity has appear at o = 2.343. Again, on the level of Hubble diagram, one 
cannot distinguish it from ACDM model. 

In this work we have also examined the significance of the Starobinsky's term. It 
appears, and can be seen on the plots of the potentials (fig 9, 11) and on the Bubble's 
diagram (fig 6), that without this term the remaining parameters of the model cannot be 
properly fitted. Quadratic gravity along well qualitatively mimics ACDM model till the 
present time (see fig. (8, 17)). However combining it with the non-minimally coupled scalar 
field provides fine tuning and much better adjustment to the concordance model. It is also to 
be observed that the parameter a is negative provided that /? / 0. In contrast for Ip=Q one 
has Q > which provides the so-called Chameleon effect [42]. We have no at the moment 
good explanation for this and the problem will be studded in our future work. 
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Figure 17. The phase portrait for quadratic gravity model //3=o on the plain (a, a) - estimations 
including CMB data (see fig. (8)). One can observe Big Bang singularity, saddle point corresponding 
to static solution as well as behaviour of the system near sudden singularity (red vertical line). Bolded 
trajectory of the flat model is determined by energy level E = 0. 



The Bayesian framework of model selection has been used for comparison theoretical 
models with the concordance ACDM one. We investigated Bayes Factor to show which 
model is the best (most probable) in the light of the astronomical data. We have found that 
while some special cases of theoretical models becomes a weak evidence in favor of it over 
ACDM if observations from current epoch are used (SN, H(z), BAO data). Inclusion of the 
information coming from early epochs (CMB data ) changes this situation because we obtain 
strong evidence in favor of Standard Cosmological Model. 

Our investigations here have aimed to distinguish the favorable model by cosmography 
of the FLRW background metric in the sample of theoretical models. Because of plenitude of 
dynamical scenarios, an introductory selection of sample of theoretical models was necessary 
and it accounts in final Bayesian inference. In the Bayesian framework adding of new obser- 
vations is natural for improving models parameters. It means that the effects of cosmological 
perturbations in this class of models have not been considered here and this important task 
is postponed for future investigations. This will allow to enlarge the discriminatory tools for 
further analysis. For example, a sound speed of the fluctuations for the quadratic gravity 
model 1/3=0 a-s calculated in [24] is 



be, in principle, rejected. The dynamics of cosmological perturbations in extended gravity 
models, including Palatini formulation, have been studied in number of papers (see e.g. 
[43] and references therein) both for a matter density as well as for the background metric. 
However the contributions coming from non-minimal dilaton-curvature couplings one deals 
here have not been developed yet (cf. [44] in the context of other scalar-tensor theories) and 
will be a subject of our future investigations. 





which in our case r^o,. 



4 yields a super luminal value > 1. Thus such a model should 
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models II: equation (5.5) - the parameters estimated without CMB data 
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models I: equation (4.6) - the parameters estimated with CMB data added 






no,/3 


,5 


Oo,™ 


90 




7 
8 
9 


a = 0, Hq^p G< -10, 10 >, (5 e (0, 1 > 

f2o,a e< 0, 40 >, ^ = 
Oo,o 6< -30, >, no,/3 e< -10, 10 >, a g (o, i > 


3.854 ± 0.103(3.843) 
-20. 824^8 815 (-6.483) 


0.0016 ± 0.0007(0.0004) 
4.237i;;?^8 (1.076) 


0.0016 ± 0.0007(0.0004) 

o.2ioi°-;;;^(o.i72) 


0.415+'^-^!'^ (0.409) 
O-awto.SiS (0-319) 
0.266 ± 0.011(0.268) 


0.355 
-1.236 
-0.844 


-0.097 
-1.157 
-0.896 


models II: equation (5.5) - the parameters estimated with CMB data added 




no, a 


no.c 


,5 




Ho,™ 


90 


»=.ff,0 


10 
11 


a = 0, no,c 6< -1, 10 >, (5 g (0, 1) 

5^0,0 6< -60, -10 >, no,c e< -1, 5 >, « e (o, i) 


-56.342+5 (-59.887) 


0.047 ± 0.005(0.046) 
1.905lJ5;227 (2.000) 


0.605^0 021 (0-604) 
0-615l:°-™3 (0.613) 


1.989to"'oi7{1.991) 
0.0121°-™! (0.012) 


0.093l'^-'^'^J (0.093) 
0.023t°-ooj (0.023) 


-0.029 
-0.453 


-0.352 
-0.635 


model ACDM: H'' / H'f^ = 1 - Oq „ + Oq „(1 + z)^ - the parameter estimated with CMB data added 






90 




12 1 no,77. e< 0, 1 > 


0.262+'„;5j2 (0.262) 


-0.608 


-0.738 



Table 3. The values of estimated parameters (mean of the marginahzed posterior probabihties and 68% credible intervals or sample square roots 
of variance, together with mode of the joined posterior probabilities, shown in brackets) for all discussed models. Model /a=o corresponds to rows 
No 1, 7; model No 2, 8; model /: No 3, 9; model IIa=o^■ No 4, 10; model //: No 5, 11; ACDM: No 6, 12. Computations were made using 

Union2 + + BAO data. We compare estimations without CMB data (top part of the table) with the one employing CMB data (bottom part). 



